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Abstract 

Numerical simulations and finite-size scaling analysis have been carried out to 
study the percolation behavior of straight rigid rods of length k (fc-mers) on two- 
dimensional square lattices. The fc-mers, containing k identical units (each one 
occupying a lattice site), were adsorbed at equilibrium on the lattice. The pro- 
cess was monitored by following the probability Ri^ki^^) that a lattice composed of 
L X L sites percolates at a concentration 9 of sites occupied by particles of size k. 
A nonmonotonic size dependence was observed for the percolation threshold, which 
decreases for small particles sizes, goes through a minimum, and finally asymptoti- 
cally converges towards a definite value for large segments. This striking behavior 
has been interpreted as a consequence of the isotropic-nematic phase transition oc- 
curring in the system for large values of k. Finally, the universality class of the 
model was found to be the same as for the random percolation model. 
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1 Introduction 



The percolation problem is a topic being increasingly considered in statistical physics. One 
reason for this current interest is that it is becoming clear that generalizations of the pure 
percolation problem are likely to have extensive applications in science and technology 
[H El El m E] • Although it is a purely geometric phenomenon, the phase transition involved 
in the process can be described in terms of an usual second-order phase transition. This 
mapping to critical phenomena made percolation a full part of the theoretical framework 
of collective phenomena and statistical physics. 

Despite of the number of contributions to this problem, there are many aspects which 
are not yet completely solved. In fact, most of the studies are devoted to the percolation of 
objects whose size coincides with the size of the lattice site (single occupancy). However, 
if some sort of correlation exists, like particles occupying several k contiguous lattice sites 
(fc-mers), the statistical problem becomes exceedingly difficult and a few studies have 
been devoted to understanding the percolation of elements occupying more than one site 

(bond) pHEiEiiiniiiiiiiaiiaiiiiiiaiisiiiz]. 

In two previous articles [131 [H], referred to as papers I and II, respectively, the per- 
colation of straight rigid rods on 2D lattices was studied. In paper I, linear fc-mers were 
deposited randomly and irreversibly. It was established that (1) the percolation thresh- 
old exhibits a exponentially decreasing function when it is plotted as a function of the 
k-mer size and (2) the problem belongs to the same universality class that the random 
percolation regardless the value of k used. Paper II was devoted to the study of aligned 
rigid rods on square lattices. In this case, the fc-mers were irreversibly deposited along 
one of the directions of the lattice. The results obtained revealed that (i) the percolation 
threshold monotonically decreases with increasing k; (ii) for any value of k [k > 1), the 
percolation threshold is higher for aligned rods than for isotropically deposited rods; (iii) 
the phase transition occurring in the system belongs to the standard random percolation 
universality class; and (iv) the intersection point of the percolation cumulant curves for 
different system sizes exhibits a nonuniversal critical behavior, varying continuously with 
changing the k-mer size. 

On the other hand, numerous experimental and numerical studies have been recently 
devoted to the analysis of equilibrium properties in systems of non-spherical particles 
[ISIIISIEDIEIIESIESIEIIESIES]. Of special interest are those studies dealing with lattice 
versions of the present problem, where the situation is much less clear than in the con- 
tinuum case. In this line, a system of straight rigid rods of length k on a square lattice, 
with two allowed orientations, was recently studied by Monte Carlo simulations [22] . The 
authors found strong numerical evidence on the existence of an isotropic-nematic (IN) 
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phase transition at intermediate densities for k > 7. The nematic phase, characterized 
by a big domain of parallel fc-mers, is separated from the isotropic phase by a contin- 
uous transition occurring at a finite density. The accurate determination of the critical 
exponents, along with the behavior of Binder cumulants, showed that the transition from 
the low-density disordered phase to the intermediate-density ordered phase belongs to the 
2D Ising universality class for square lattices and the three-state Potts universality class 
for honeycomb and triangular lattices [231 |2l]. addition, the comparison between the 
configurational entropy of the system and that corresponding to a fully aligned system 
[25] (i) confirmed the existence of a IN phase transition at intermediate densities for 
long rods, (ii) allowed us to estimate the minimum value of k, which leads to the formation 
of a nematic phase, and (iii) provided an interesting interpretation of this critical value. 

Later, the dependence of the critical density on the magnitude of the lateral interac- 
tions was studied for a system of attractive rigid rods on square lattices with two allowed 
orientations EB] ■ The obtained results revealed that the orientational order survives 
in a wide range of lateral interactions. 

A notable feature is that nematic order is only stable for sufficiently large aspect 
ratios while systems of short rods do not show nematic order at all. The long-range 
orientational order also disappears in the case of irreversible adsorption (as is the case 
of papers I and II), where the distribution of adsorbed objects is different from that 
obtained at equilibrium. Thus, the irreversible adsorption leads to intermediate states 
characterized by an isotropic distribution of the directions of the adsorbed molecules, and 
to a final state (known as jamming state), in which no more objects can be deposited due 
to the absence of free space of appropriate size and shape (the jamming state has infinite 
memory of the process and the orientational order is purely local) |29j . 

It is of interest then to study the percolation problem of /c-mers at equilibrium, not only 
because of this problem has not been treated in the literature, but also because, depending 
on the values of the density and the size k, the model leads either to an isotropic system 
(as studied in paper I) or an anisotropic system (as studied in paper II). In this sense, 
the aim of the present work is to study, by Monte Carlo (MC) simulations and finite-size 
scaling analysis, the effect of the orientational order on the percolation properties for a 
system of rigid rods adsorbed at equilibrium on two-dimensional square lattices. 

2 Model and Monte Carlo simulation 

The calculations have been developed for straight rigid rods of length k, with k ranging 
between 2 and 12. Small adsorbates with spherical symmetry would correspond to the 
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monomer limit {k = 1). The distance between k-mer units is assumed to be equal to 
the lattice constant; hence exactly k sites are occupied by a k-mer when adsorbed. The 
surface was represented as a 2D square lattice of M = Lx L adsorptive sites with periodic 
boundary conditions. The only interaction between different rods is hard-core exclusion: 
no site can be occupied by more than one fc-mer. 

To describe a system of fc-mers adsorbed on M sites at a given temperature T 
and chemical potential /i, the occupation variable cj was introduced {cj = or 1, if the 
site j is empty or occupied by a A;-mer unit, respectively). Then the adsorbed phase is 
characterized by the Hamiltonian H = (cq ~ f^)Ylj'^j^ where the sum run over the M 
sites and eo is the adsorption energy of a k-mei unit (in the simulations, eo is set equal to 
zero without any loss of generality). 

The thermodynamic equilibrium is reached in the grand canonical ensemble by using 
the hyper-parallel tempering Monte Carlo (HPTMC) simulation method [SUIEI] combined 
with an efficient cluster algorithm (the pocket algorithm) [S2]. 

The HPTMC method consists in generating a compound system of R noninteracting 
replicas of the system under study. The i-th replica is associated with a chemical potential 
/ij. To determine the set of chemical potentials, {/ij}, we set the lowest chemical potential, 
/ii, in the isotropic phase where relaxation (correlation) time is expected to be very 
short and there exists only one minimum in the free energy space. On the other hand, 
the highest chemical potential, fiR, is set in the nematic phase whose properties we are 
interested in. Finally, the difference between two consecutive chemical potentials, /ij 
and /ij+i with yUj < fii+i, is set as Afi = {fii — fiR) /(i? — 1) (equally spaced chemical 
potentials). The parameters used in the present study were as follows: R = 40, fii = —2.5 
and fiR = 2.375. With these values of the chemical potential, the corresponding values 
of the surface coverage varied from 6*1 (/ii) ~ 0.21 to Or^hr) ~ 0.88 for k = 2, and from 
^i(/xi) ^ 0.33 to enifXR) ^ 0.81 for k = 12. 

Under these conditions, the algorithm to carry out the simulation process is built on 
the basis of two major subroutines: replica-update and replica- exchange. 
Replica-update: The adsorption-desorption and pocket procedure is as follows: (1) One 
out of R replicas is randomly selected. (2) A linear fc-uple of nearest-neighbor sites, 
belonging to the replica selected in (1), is chosen at random. Then, if the k sites are 
empty, an attempt is made to deposit a rod with probability W = min{l, exp(/3/i)}; if 
the k sites are occupied by units belonging to the same k-mer, an attempt is made to 
desorb this k-mer with probability W = min{l, exp(— /3/i)}; and otherwise, the attempt 
is rejected. In addition, the pocket algorithm [32] is applied in order to reach equilibrium 
in a reasonable time. This algorithm includes the following steps: (1) Choose randomly 
a lattice symmetry axis and a seed fc-mer. The seed k-mer now is the sole element of a 
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set V and the rest of fc-mers belongs to a set O. (2) Move elements from the set V to O 
in the following way: Pick up an arbitrary element i oiV and reflect with respect to the 
symmetry axis. If i overlaps with other elements of (9, the latter are transferred from O 
to V. (3) Repeat the step (2) until the set V = 0. 

Replica- exchange: Exchange of two configurations Xj and Xj, corresponding to the i-th 
and j-th replicas, respectively, is tried and accepted with probability W = min{l, exp (—A) } 
Where A in a nonthermal grand canonical ensemble is given by [—f3{fij — fii) {Nj — Ni)], 
and Ni (Nj) represents the number of particles of the i-th (j-th) replica. 

The complete simulation procedure is the following: (1) replica- up date, (2) replica- 
exchange, and (3) repeat from step (1) RM times. This is the elementary step in the 
simulation process or Monte Carlo step (MCs). 

For each value of the chemical potential fii, the equilibrium state can be well repro- 
duced after discarding the first tq MCs. Then, a set of m samples in thermal equilibrium is 
generated. The corresponding surface coverage 9i{fii) is obtained through simple averages 
over the m samples (m MCs). 

As mentioned before in Ref. [22], the relaxation time increases very quickly as the 
A;-mer size increases. Consequently, MC simulations for large adsorbates are very time 
consuming and may produce artifacts related to non-accurate equilibrium states. In order 
to discard this possibility, equilibration times tq of the order 0(10^ MCs) were used in 
this study, with an effort reaching almost the limits of our computational capabilities [33]. 




3 Finite-size scaling results 

The central idea of the pure percolation theory is based in finding the minimum con- 
centration of elements (sites or bonds) for which a cluster extends from one side to the 
opposite one of the system. At this particular value of the concentration, the percolation 
threshold 6c, a second-order phase transition occurs in the system, which is characterized 
by well-defined critical exponents [2]. 

The finite-size scaling theory [34J gives us the basis to achieve the percolation threshold 
and the critical exponents of a system with a reasonable accuracy. For this purpose, the 
probability R = Rf^f^{9) that a L x L lattice percolates at a concentration p of sites 
occupied by rods of size k can be defined [21 [35]. Here, the following definitions can 
be given according to the meaning of X: a) R^^^\6) = the probability of finding a 
rightward (downward) percolating cluster; and b) R^j.{6) = the probability of finding 
either a rightward or a downward percolating cluster. 

In the Monte Carlo (MC) simulations, each run consists of the following steps: (a) 
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the construction of m = 2 x 10^ samples for the desired fraction 6{fi) of occupied sites 
(according to the HPTMC method described above); and (b) the cluster analysis by using 
the Hoshen and Kopelman algorithm [36]. In the last step, the existence of a percolating 
island is verified for each sample. The spanning cluster could be determined by using the 
criteria R, D and U . This scheme is carried out for obtaining the number of samples 
for which a percolating cluster of the desired criterion X is found. Then, R^i.{0) = /m 
is defined and the procedure is repeated for different values of 9, L and k. 

Besides the parameter R, the percolation order parameter P = {Sl)/L'^ [371 |38] has 
been measured, where Sl represents the size of the largest cluster and (...) means the 
average over the m samples. The corresponding percolation susceptibility x has also been 
calculated, x = [{SD - (Sl)^] /L^. 

In Figure 1, the probability Rl ^{9) is presented for /c = 3 (left), k = 5 (middle) 
and k = 7 (right) and different lattice sizes L/k = 5,10,15,20 [32j, as indicated in the 
caption of the figure. Several conclusions can be drawn from the figure. First, curves cross 
each other in a unique point R^ (measured in the vertical axes), which depends on the 
criterion X used (data corresponding to criteria R and D are not shown in the figure for 
clarity) and those points are located at very well defined values in the 6'-axes determining 
the critical percolation threshold (measured in the horizontal axes) for each k. The values 
of 9c{k) obtained by following the procedure of Figure 1 are collected in Table 1 and are 
plotted in Figure 2a (full circles). A nonmonotonic size dependence is observed for the 
percolation threshold, which decreases for small particles sizes, goes through a minimum 
around k = 5, and finally asjTiiptotically converges towards a definite value for large 
segments. This striking behavior can be interpreted as a consequence of the IN phase 
transition occurring in the system for large values of k [221 ESI EH |25l |26]. To understand 
this effect, it is convenient first to recall some results about percolation of straight rigid 
rods isotropically deposited (or isotropic k-mers) [I3j and straight rigid rods deposited 
along one of the directions of the lattice (or aligned fc-mers) [T7] . 

In the first case (isotropic fc-mers) |13j, the rods are deposited randomly and irre- 
versibly. The process is as follows: a linear fc-uple of nearest-neighbor sites is randomly 
selected; if it is vacant, a fc-mer is then deposited on those sites. Otherwise, the attempt 
is rejected. This filling process, known as Random Sequential Adsorption model, leads to 
states characterized by an isotropic distribution of the directions of the adsorbed fc-mers. 
In the case of aligned fc-mers fUl, the adsorption process is also irreversible, but this time 
the k-mers are deposited along one of the directions of the lattice, forming a nematic 
phase. 

As it can be observed from Figure 2a, the percolation threshold curves corresponding 
to isotropic (open circles) and aligned (open squares) fc-mers decrease upon increasing k. 
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At the beginning, for small values of k, the curves rapidly decrease. However, they flatten 
out for larger values of k and finally asymptotically converge towards a definite value 6* 
(as k — 7- oo), being 6* ~ 0.46 for isotropic fc-mers and 6* ~ 0.54 for aligned fc-mers. 

Now, the nonmonotonic behavior observed in Figure 2a can be explained on the basis 
of the results for isotropic and aligned fc-mers. In fact, for small values of k, there is 
no orientational order in the equilibrium adlayer and Oc{k) decreases following the curve 
corresponding to isotropic fc-mers. On the other hand, for large values of k, a nematic 
phase, characterized by a big domain of parallel fc-mers, is formed on the surface at an 
intermediate density Ojn [221 [221 1211 12S1 [2S]- The values of Ojiy reported in Ref. are 
listed in the third column of Table 1. For k > 9, 6^ > 9 in and, consequently, (1) the 
percolation properties of the orientational phase are similar to the ones corresponding to 
perfectly aligned fc-mers and (2) 9c{k) follows the curve corresponding to aligned /c-mers. 
As a consequence of the behavior of Oc{k) in the limits of small and large fc-mers, a marked 
minimum is observed at intermediate values of k. 

Typical equilibrium states are shown schematically in parts (b)-(d) of Figure 2. In 
the case of Figure 2b [point i in Figure 2a: k = 2 and 6 ~ 0.53], the equilibrium state 
corresponds to an isotropic and nonpercolating state. The situation is different in Figure 
2c [point / in Figure 2a: k = 5 and 6 ^ 0.53], where the equilibrium state is an isotropic 
and percolating state. Thin bars represent isolated /c-mers and thick bars (red bars in 
the online version) correspond to fc-mers belonging to the percolating cluster. Finally, a 
nematic and nonpercolating state is observed in Figure 2d [point j in Figure 2a: k = 10 
and 9 ^ 0.53]. 

A second conclusion of Figure 1 is that the crossing points of the percolation cumulants 
i?^* (as Rf and Rf*) vary as the k-mei size is increased. This behavior was already 
observed in the case of perfectly aligned fc-mers irreversibly deposited on 2D lattices [T^ 
and can be attributed to the strong anisotropy introduced in the system by the IN phase 
transition occurring for large fc-mers. In the case of thermal transitions, the link between 
anisotropy and a nonuniversal behavior of the intersection point of the cumulants has 
been discussed in [151 US]- In fact, as pointed out by Selke et al. [13, SH], the measure of 
the cumulant intersection may depend on various details of the model, which do not affect 
the universality class, in particular, the boundary condition, the shape of the lattice, and 
the anisotropy of the system. In order to confirm or discard this hypothesis, the critical 
exponents u, (i and 7 were calculated as follows. 

The value of v can be obtained through the scaling relationship for R^f,{9) 



(1) 



being R^{u) the scaling function and u 



{6 — 6c) L^^'^ . Then, the maximum of the 
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derivative of Eq. leads to {dR^,^/de)^^^ oc L^/''. 

In the inset of Figure 3, [dR^j^/dO)^^^ has been plotted as a function of L/k (note 
the log-log scale) for different fc-mers as indicated. According to Eq. ([1]) the slope of each 
line corresponds to As it can be observed, the slopes of the curves remain constant 
(and close to 3/4) for all values of k. The results coincide, within numerical errors, with 
the exact value of the critical exponent of the ordinary percolation v = 4/3. 

The scaling behavior can be further tested by plotting R^^.{0) vs {6 — 9c) L^^" and 
looking for data collapsing. Using the values of 6c previously calculated and the exact 
value 1/ = 4/3, an excellent scaling collapse was obtained (Figure 3) for all value of k-mer 
size. This leads to independent control and consistency check of the numerical value of 
the critical exponent u. 

The critical exponents /3 and 7 were obtained from the scaling behavior of P and x 

p = L-^^'^p [\e - ec\L^'''] , (2) 

and 

X = L-^/^X [{0 - 0c) L''"] , (3) 

where P and x are scaling functions for the respective quantities. According to Eqs. 
02]) and ([3]), Figure 4 shows the collapse of the curves of P and x (inset) for a typical 
/c-mer size {k = 5) and different lattice sizes as indicated. The study was repeated for all 
studied values of k (these data are not shown here for lack of space). The data scaled 
extremely well using the exact percolation exponents /3 = 5/36 and 7 = 43/18. The 
results obtained in Figures 3 and 4 clearly supports that, as in the case of isotropic rods 
|13] . the universality class of the problem does not depend on the k-mer size. 

4 Conclusions 

In this work, Monte Carlo simulations and finite-size scaling theory have been used to 
study the percolation properties of straight rigid rods of length k adsorbed at equilibrium 
on two-dimensional square lattices. 

A nonmonotonic size dependence was found for the percolation threshold 6c, which 
decreases for small particles sizes [following the curve of 6c{k) for /c-mers irreversibly and 
isotropically deposited], goes through a minimum, and finally asymptotically converges 
towards a definite value for large segments [following the curve of 6c{k) for aligned /c-mers 
irreversibly deposited along one of the directions of the lattice]. This striking behavior 
has been interpreted as a consequence of the IN phase transition occurring in the system 
for large values of k. Finally, the analysis of the critical exponents u, /3 and 7 revealed 
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that the percolation phase transition involved in the problem considered in the present 
paper belongs to the same universality class of the ordinary random percolation. 
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Table and Figure Captions 



Table 1: Values of the percolation threshold 6c for k ranging from 1 to 12 and of the 
isotropic-nematic critical density 6 in [25] for k between 7 to 12. Note that the IN phase 
transition on a square lattice with two allowed orientations occurs only for k > 7 [221 
l23t |2H |25| |26] . In the case of = 1 , the problem of percolation of monomers on square 
lattices has been one of the most studied percolation models in the literature, and the 
percolation threshold has been measured multiple times as 0.592746a;x, where the last 
two decimal places (and error) are 21(13) [40], 21(33) [HJ, 03(09) |2], and 06(05) 
For a complete overview on this topic we refer the reader to Ref. [H]. 

Fig. 1: Fraction of percolating lattices R^f. as a function of the concentration 6' for A; = 3 
(left), k = 5 (middle) and k = 7 (right) and different lattice sizes: L/k = 5, diamonds; 
L/k = 10, squares; L/k = 15, circles and L/k = 20, triangles. Horizontal dashed lines 
show the intersection points R^* ■ Vertical dashed line denotes the percolation threshold 
6c in the thermodynamic limit. 

Fig. 2: (Color online) (a) The percolation threshold 6c as a function of k for straight rigid 
rods on 2D square lattices and three different adsorption processes: nematic irreversible 
adsorption (open squares), isotropic irreversible adsorption (open circles) and equilibrium 
adsorption (solid circles). The dashed lines are simply a guide for the eye. In all cases, 
the error bar is smaller that the size of the symbols. The asymptotic limits, 6*, [6c{k) 
for — 7- oo] are shown, (b) Schematic representation of the equilibrium state for k = 2 
and 6 ~ 0.53 (point i in the figure), (c) Same as (b) for k = 5 and 6 ^ 0.53 (point / 
in the figure). Thin bars represent isolated /c-mers and thick bars (red bars in the online 
version) correspond to /c-mers belonging to the percolating cluster, (d) Same as (b) for 
A; = 10 and 6 ~ 0.53 (point j in the figure). 

Fig. 3: Data collapsing of the fraction of percolating samples Rli^{6) as a function 
of the argument {6 — 6c) L^^'^. Each set of curves corresponds to a different value of k 
as indicated. For each k, different lattice sizes {L/k = 5, 10, 15 and 20) have been 
considered. Inset: In ( — j as a function of ln(L/A;) for different values of k as 

V / max 

indicated. According to Eq. ([1]) the slope of each line corresponds to l/u. 

Fig. 4: Data collapsing of the order parameter, PL^^'^ vs \6 — 6c\L^^'^, and of the suscep- 
tibility, xL'^'^" vs {6 - 6c) L'^'" (inset) for = 5. The plots were made using the exact 
percolation exponents z/ = 4/3,/3 = 5/36 and 7 = 43/18. 
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Table 1 



7 • 7 

k — mer size, k 


1 j_' j_l 111/1 

percolation threshold, 


isotropic — nematic critical density, c^/jv 


1 


0.592746xx(yy) 




2 


0.586(1) 




3 


0.575(1) 




4 


0.567(1) 




5 


0.563(1) 




6 


0.559(1) 




7 


0.593(1) 


0.729(3) 


8 


0.586(1) 


0.648(3) 


9 


0.575(1) 


0.569(3) 


10 


0.567(1) 


0.502(3) 


11 


0.563(1) 


0.457(5) 


12 


0.559(1) 


0.413(5) 
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Fig. 1: D. A. Matoz-Fernandez et al . 
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Fig. 2a: D. A. Matoz-Fernandez et al . 



16 



Point i,k = 2 



i- i ii' zz 

— 1 1 1 _ — — < 
= 1 I _I I _ I I _ _ I I - _ I 

ri'i "i= 'iJ-"- " lE 
1 r'T-"' i" -I ' "i"' 

T ' M ' I i' T - = 

_ ' I ri.iii I - - I ~ 

— |i __ i -| - II I 

=!_ v"iiip""i' ' r ! ,--^=11,' 

- =r -- ' I r- m' '- 
i-j":l=_i_ ! iz-r'i" 

_-_-i-_ I _ii' i-_ji-_ II I- 

I ""i' " I r' I ' " 

'ii'i — ' ' — 'I'l I 

MINI ~l ~l ~ ~ 



(b) 



Fig. 2b: D. A. Matoz-Fernandez et al . 
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Fig. 2c: D. A. Matoz-Fernandez et al . 
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Point j,k = 10 
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Fig. 2d: D. 



Matoz-Fernandez et al . 
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Fig. 3: D. A. Matoz-Fernandez et al . 
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Fig. 4: D. A. Matoz-Fernandez et al. 
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